set linesize 140
set rmsg on
set scheme s1mono
/***************************************************************************
	Project		:	Army Service in the All-Volunteer Era	
	Author(s)	:	Kyle Greenberg	(kyle.greenberg@westpoint.edu)
					Matthew Gudgeon (matthew.gudgeon@westpoint.edu)
					Adam Isen 		(Adam.Isen@treasury.gov)
					Corbin Miller 	(Corbin.Miller@treasury.gov)
					Rich Patterson 	(rich_patterson@byu.edu)
	File Name	:	robust.do
	Description	:	This file runs robustness checks for main analysis to 
					bandwidths, controls, and functional form.
****************************************************************************/

*---- if not running master, set up file structure ----*
*change directory to where programs and subfolders are stored
cd ""

if "${raw}"==""		global raw		"raw/"
if "${data}"==""	global data		"data/"
if "${output}"==""	global output	"output/"

cap mkdir	"${output}"
cap mkdir	"${output}robust/"

*---- select subprograms ----*
local robust_regs 1
local robust_figs 1

*---------------------*
*---- regressions ----*
*---------------------*
if `robust_regs'==1{
	local outcome wages
	foreach strata in all black white {
		*initialize output file
		clear
		set obs 17
		gen bw = _n + 2
		foreach cont in nocont cont {
			foreach model in lin lintri quad quadtri {
				foreach cut in 31 50 {
					foreach x in m b se N{
						qui gen `x'_`cut'_`model'_`cont' = .
					}
				}
			}
		}
		save "${output}robust/robust_`strata'", replace

		forval bw=3/19{
			foreach cut in 31 50 {
				foreach control in nocont cont{
					use "${data}cumulative-analysis" if `strata'==1, clear
					gen wt`cut'_`bw' = 1-abs((firstafqt-`cut'+.5)/`bw')
					if "`control'" == "nocont"{
						local controls i.quarter
						local controllist quarter
					}
					if "`control'" == "cont"{
						local controls i.quarter i.app_state age_days male black hisp educ_lessthanhs educ_inhs educ_ged educ_some_coll educ_coll_grad
						local controllist quarter app_state age_days male black hisp  educ_lessthanhs educ_inhs educ_ged educ_some_coll educ_coll_grad
					}

					keep `outcome'_11_19_mean inst`cut' k`cut' k`cut'_2 instk`cut' instk`cut'_2 `controllist' access wt`cut'_`bw' firstyear firstafqt
					
					keep if inrange(firstafqt,`cut'-`bw',`cut'+`bw'-1)

					*Linear (No Kernel)
					ivregress 2sls `outcome'_11_19_mean k`cut' instk`cut' `controls' (access=inst`cut'), robust
					local b_`cut'_lin_`control' = _b[access]
					local se_`cut'_lin_`control' = _se[access]
					local N_`cut'_lin_`control' = e(N)
					sum `outcome'_11_19_mean if e(sample)
					local m_`cut'_lin_`control' = r(mean)
					
					*Linear (Triangle Kernel)
					ivregress 2sls `outcome'_11_19_mean k`cut' instk`cut' `controls' (access=inst`cut') [aw=wt`cut'_`bw'], robust
					local b_`cut'_lintri_`control' = _b[access]
					local se_`cut'_lintri_`control' = _se[access]
					local N_`cut'_lintri_`control' = e(N)
					sum `outcome'_11_19_mean if e(sample)
					local m_`cut'_lintri_`control' = r(mean)
					
					*Quadratic (No Kernel)
					ivregress 2sls `outcome'_11_19_mean k`cut' k`cut'_2 instk`cut' instk`cut'_2 `controls' (access=inst`cut'), robust
					local b_`cut'_quad_`control' = _b[access]
					local se_`cut'_quad_`control' = _se[access]
					local N_`cut'_quad_`control' = e(N)
					sum `outcome'_11_19_mean if e(sample)
					local m_`cut'_quad_`control' = r(mean)
					
					*Quadratic (Triangle Kernel)
					ivregress 2sls `outcome'_11_19_mean k`cut' k`cut'_2 instk`cut' instk`cut'_2 `controls' (access=inst`cut') [aw=wt`cut'_`bw'], robust
					local b_`cut'_quadtri_`control' = _b[access]
					local se_`cut'_quadtri_`control' = _se[access]
					local N_`cut'_quadtri_`control' = e(N)
					sum `outcome'_11_19_mean if e(sample)
					local m_`cut'_quadtri_`control' = r(mean)
					
					*STORE COEFFICIENT ESTIMATES 
					use "${output}robust/robust_`strata'", clear
					foreach model in lin lintri quad quadtri{
						foreach x in b se N m{
							di "`x'_`cut'_`model'_`control' = ``x'_`cut'_`model'_`control''"
							replace `x'_`cut'_`model'_`control' = ``x'_`cut'_`model'_`control'' if bw == `bw'
						}
					}
					save  "${output}robust/robust_`strata'", replace
				}
			}
		}
	}
}	//end if
else di "*---- skipping robustness regressions ----*"

if `robust_figs'==1{
	local outcome wages
	foreach strata in all black white {
		use "${output}robust/robust_`strata'", clear
		gen y1 = bw-.1
		gen y2 = bw+.1		
		foreach cut in 31 50 {
			foreach model in lin lintri quad quadtri {
				foreach control in cont nocont {
					gen hi_`cut'_`model'_`control' = b_`cut'_`model'_`control' + 1.96*se_`cut'_`model'_`control'
					gen lo_`cut'_`model'_`control' = b_`cut'_`model'_`control' - 1.96*se_`cut'_`model'_`control'
				}
				tw 	(rcap hi_`cut'_`model'_nocont lo_`cut'_`model'_nocont y1, lcolor(gs7) lp(longdash) lw(medthick))	///
					(scatter b_`cut'_`model'_nocont y1, mcolor(gs7) msize(vsmall))	///		
					(rcap hi_`cut'_`model'_cont lo_`cut'_`model'_cont y2, lcolor(black) lp(shortdash) lw(medthick))	///
					(scatter b_`cut'_`model'_cont y2, mcolor(black) msize(vsmall))	///				
					,  yline(0, lpattern(shortdash))	///
					legend( order(2 "Baseline" 4 "Added Controls") )	///
					ytitle("Earnings (Dollars)") ylabel(,angle(0) format(%9.0fc))	///
					xtitle("Bandwidth") xtick(3(1)19) xlabel(3(1)19)
				graph export "${output}robust/robust_to_bw_`strata'_`cut'_`model'.pdf", as(pdf) replace
			}
		}
	}
}	//end if
else di "*---- skipping robustness figures ----*"
